function hv = Hvec(invP,etaModel,j,beta,WmtWm,v)

Jv = Jvec_TDIP(invP,etaModel,j,v);
JtJv = JTvec_TDIP(invP,etaModel,j,Jv);

hv = JtJv + beta*WmtWm*v;
